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Abstract.  In  this  paper,  we  study  the  strong  stability  preserving  (SSP)  prop¬ 
erty  of  a  class  of  deferred  correction  time  discretization  methods,  for  solving  the 
method-of-lines  schemes  approximating  hyperbolic  partial  differential  equations. 


1.  Introduction 

In  this  paper,  we  are  interested  in  the  numerical  solutions  of  hyperbolic  partial 
differential  equations  (PDEs).  A  typical  example  is  the  nonlinear  conservation  law 

(1.1)  ut  =  ~f(u)x. 

A  commonly  used  approach  to  design  numerical  schemes  for  approximating  such 
PDEs  is  to  first  design  a  stable  spatial  discretization,  obtaining  the  following  method- 
of-lines  ordinary  differential  equation  (ODE)  system 

(1.2)  ut  =  L(u ) 

to  approximate  (1.1).  Notice  that  even  though  we  use  the  same  letter  u  in  (1.1)  and 

(1.2) ,  they  have  different  meanings.  In  (1.1),  u  =  u(x,t)  is  a  function  of  x  and  t, 
while  in  (1.2),  u  =  u(t )  is  a  (vector)  function  of  t  only.  Stable  spatial  discretization 
for  (1.1)  includes,  for  example,  the  total  variation  diminishing  (TVD)  methods  [6], 
the  weighted  essentially  non-oscillatory  (WENO)  methods  [7],  and  the  discontinuous 
Galerkin  (DG)  methods  [1].  In  this  paper,  we  assume  that  the  spatial  discretization 

(1.2)  is  stable  for  a  first  order  Euler  forward  time  discretization 

(1.3)  un+l  =  un  +  A  tL(un) 
under  a  suitable  time  step  restriction 

(1.4)  At  <  A t0. 

This  stability  is  given  as 

(1.5)  |K+1||  <  |K|| 
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for  a  suitable  norm  or  semi-norm  ||  •  ||.  For  the  TVD  schemes  [6],  ||  •  ||  is  taken  as 
the  total  variation  semi-norm.  For  technical  reasons,  we  would  also  need  a  different 
but  closely  related  spatial  discretization  to  (1.1): 

(1.6)  ut  =  L(u) 

with  the  property  that  the  following  first  order  “backward”  time  discretization 

(1.7)  un+l  =  un-  A tL(un) 

is  stable  in  the  sense  of  (1.5)  under  the  same  time  step  restriction  (1.4).  For  the 
conservation  law  (1.1),  the  operator  L  can  often  be  obtained  simply  by  reversing 
the  wind  direction  in  the  upwind  approximation.  We  refer  to,  e.g.  [11],  [7]  and  [1] 
for  such  implementation  in  ENO,  WENO  and  DG  methods. 

Even  though  the  fully  discretized  scheme  (1.3)  is  assumed  to  be  stable  as  in  (1.5), 
it  is  only  first  order  accurate  in  time.  For  a  high  order  spatial  discretization  such 
as  in  the  WENO  and  DG  methods,  we  would  certainly  hope  to  have  higher  order 
accuracy  in  time  as  well.  A  higher  order  time  discretization  for  (1.2)  is  called  strong 
stability  preserving  (SSP)  with  a  CFL  coefficient  c,  if  it  is  stable  in  the  sense  of  (1.5) 
under  a  possibly  modified  time  step  restriction 

(1.8)  At  <  cAt0. 

SSP  time  discretizations  were  first  developed  in  [10]  for  multi-step  methods  and  in 
[11]  for  Runge-Kutta  methods.  They  were  referred  to  as  TVD  time  discretizations 
in  these  papers,  since  the  semi-norm  involved  in  the  stability  (1.5)  was  the  total 
variation  semi-norm.  More  general  SSP  time  discretizations  can  be  found  in,  e.g. 
[4,  12,  13,  3].  The  review  paper  [5]  summarizes  the  development  of  the  SSP  method 
until  the  time  of  its  publication. 

In  this  paper  we  study  the  SSP  property  of  a  newly  developed  time  discretization 
technique,  namely  the  (spectral)  deferred  correction  (DC)  method  constructed  in 
[2],  An  advantage  of  this  method  is  that  it  is  a  one  step  method  (namely,  to  march 
to  time  level  n  +  1  one  would  only  need  to  store  the  value  of  the  solution  at  time 
level  n)  and  can  be  constructed  easily  and  systematically  for  any  order  of  accuracy. 
This  is  in  contrast  to  Runge-Kutta  methods  which  are  more  difficult  to  construct 
for  higher  order  of  accuracy,  and  to  multi-step  methods  which  need  more  storage 
space  and  are  more  difficult  to  restart  with  a  different  choice  of  the  time  step  At. 
Linear  stability,  such  as  the  A-st  ability,  A(«)-stability,  or  L-stability  issues  for  the 
DC  methods  were  studied  in,  e.g.  [2,  8,  14].  However,  for  approximating  hyperbolic 
equations  such  as  (1.1)  with  discontinuous  solutions,  linear  stability  may  not  be 
enough  and  one  would  hope  the  time  discretization  to  have  the  SSP  property  as 
well. 

The  (s  +  l)-th  order  DC  time  discretization  to  (1.2)  that  we  consider  in  this 
paper  can  be  formulated  as  follows.  We  first  divide  the  time  step  [tn,tn+1]  where 
tn+i  =  into  s  subintervals  by  choosing  the  points  tJm^  for  m  —  0, 1,  •  •  •  ,  s  such 

that  tn  =  <  t W  <  •  •  •  <  f(m)  <  . . .  <  f(s)  =  tn+1.  We  use  Af(m)  =  £(m+1)  —  £(m)  to 

denote  the  sub-time  step  and  u^")  to  denote  the  k-th  order  approximation  to  u(t ^). 
The  nodes  Am)  can  pe  chosen  equally  spaced,  or  as  the  Chebyshev  Gauss-Lobatto 
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nodes  on  [tn,tn+1]  for  high  order  accurate  DC  schemes  to  avoid  possible  instability 
associated  with  interpolation  on  equally  spaced  points.  Starting  from  un,  the  DC 
algorithm  to  calculate  un+l  is  in  the  following. 

Compute  the  initial  approximation: 

(0)  n 
U[  J  =  U  . 

Use  the  forward  Euler  method  to  compute  a  first  order  accurate  approxi¬ 
mate  solution  u \  at  the  nodes  {t^}^l=1: 

For  m  —  0,  •  •  •  ,  s  —  1 

(1.9)  4m+1)  =  +  At^Liu^). 


Compute  successive  corrections: 

For  k  —  1,  •  •  •  ,  s 

(0)  n 
Uk+ 1  =  U  ■ 

For  m  —  0,  •  •  •  ,  s  —  1 

(1.10)  4« 11  =  45  +  »tAi<"»(i(4”> )  -  L(4"*>))  +  C +1(L(ut)), 

where 

(1.11)  0  <  6k  <  1 

and  /™+1(L(tifc))  is  the  integral  of  the  s-th  degree  interpolating  polynomial  on  the 
s  +  1  points  (t^\  L('ujy>))|=0  over  the  subinterval  [t^m\  which  is  the  numerical 

quadrature  approximation  of 


(1.12) 

Finally  we  have  un+1  = 


L{u[r))dT. 


The  scheme  described  above  with  6 *.  =  1  is  the  one  discussed  in  [2,  8].  In  [14],  the 
scheme  is  also  discussed  with  general  0  <  9k  <  1  to  enhance  linear  stability.  The 
term  with  the  coefficient  6k  does  not  affect  accuracy. 

In  the  next  three  sections  we  will  study  the  SSP  properties  of  the  DC  time  dis¬ 
cretization  for  the  second,  third  and  fourth  order  accuracy  (s  =  1,  2,  3),  respectively. 
In  section  5  we  will  provide  a  numerical  example  of  using  the  SSP  DC  time  dis¬ 
cretizations  coupled  with  a  WENO  spatial  discretization  [7]  to  solve  the  Burgers 
equation.  Concluding  remarks  are  given  in  section  6. 


2.  Second  order  discretization 

For  the  second  order  (s  =  1)  DC  time  discretization,  there  is  no  subgrid  point 
inside  the  interval  \tn)tn+1].  We  can  easily  work  out  the  explicit  form  of  the  scheme 

41}  =  un  +  A  tL(un) 
un+ 1  =  un  +  t  ( L(un )  +  L(ui1})) 


(2.1) 
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Notice  that  this  is  exactly  the  optimal  second  order  SSP  Runge-Kutta  scheme  orig¬ 
inally  given  in  [11]  and  proven  optimal  for  the  SSP  property  among  all  second  order 
Runge-Kutta  schemes  in  [4].  The  CFL  coefficient  c  in  (1.8)  for  this  scheme  is  1. 

Even  though  the  SSP  property  for  the  scheme  (2.1)  was  already  proven  in  [11,  4], 
we  will  prove  it  again  here  to  illustrate  the  approach  that  we  will  use  also  for  higher 
order  DC  time  discretizations.  This  approach  was  used  in  [12]  to  study  SSP  Runge- 
Kutta  methods.  The  first  equation  in  (2.1)  is  already  in  Euler  forward  format.  The 
idea  of  the  proof  is  to  write  the  second  equation  in  (2.1)  as  a  convex  combination  of 
Enlcr  forward  steps.  That  is,  for  arbitrary  oq,  a2  satisfying 

(2.2)  T  0,  a2  T  0,  a i  T  a2  =  1, 

we  rewrite  the  second  equation  in  (2.1)  as 

un+1  =  alUn  +  ^A  tL(un)  +  a2un  +  ^A  tL(u^]) 

and  substitute  the  first  equation  in  (2.1)  into  the  a2un  term  of  the  equation  above 
to  obtain 

(1  2a  \  /  ^ 

un  H - - — —A tL{un)  )  +  a2  \  u<L1)  +  - — A tL(u^) 

2a  1  J  \  2a  2 

Clearly,  this  is  a  convex  combination  of  two  Euler  forward  steps.  By  assumption, 
the  first  order  Euler  forward  step  (1.3)  is  stable  in  the  sense  of  (1.5)  under  the  time 
step  restriction  (1.4),  hence  it  is  clear  that  (2.3)  is  stable  in  the  sense  of  (1.5)  under 
the  modified  time  step  restriction 

1  -  2c*2  1 

- -At  <  A  to,  - At  <  A  t0. 

2a  i  ~  ’  2a2  “ 

Notice  that  aq,  a2  are  arbitrary  subject  to  (2.2),  hence  the  CFL  coefficient  c  defined 
in  (1.8)  for  the  step  (2.3),  hence  the  scheme  (2.1),  to  be  SSP  is 

f  2a  i 

(2.4)  c  =  max  min  <  - ,  2a2 

(  1  —  2  a2 

where  the  optimization  is  taken  subject  to  the  constraint  (2.2).  As  in  [12],  we 
reformulate  the  optimization  problem  (2.4)  as 

(2.5)  c  =  max  z 

{01,02} 

subject  to  the  constraint  (2.2)  and 

(2.6)  2a  1  >  z(  1  —  2 a2),  2a2  >  z. 

We  then  use  the  Matlab  routine  “fminicon”  to  obtain  the  solution  c.  The  Matlab 
routine  produces  the  optimal  solution  c  =  1  achieved  at  a\  —  a2  —  \.  This  is  the 
same  result  as  the  one  already  obtained  in  [4]  theoretically.  Of  course,  for  this  simple 
optimization  problem,  it  is  not  necessary  to  use  the  Matlab  routine.  However  for 
the  more  complicated  optimization  problems  later  associated  with  higher  order  DC 
schemes,  the  usage  of  this  Matlab  routine  will  be  helpful. 
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We  remark  that  the  sole  purpose  of  writing  the  second  equation  of  (2.1)  into  the 
mathematically  equivalent  but  more  complicated  form  (2.3)  is  to  obtain  the  optimal 
CFL  coefficient  c  in  (1.8)  for  the  provable  SSP  property  of  the  scheme  (2.1).  In 
actual  computation  we  would  use  (2.1)  since  it  is  simpler  to  implement. 


3.  Third  order  discretization 


For  the  third  order  (s  =  2)  DC  time  discretization,  there  is  only  one  subgrid 
point  inside  the  interval  [tn,tn+1].  By  symmetry,  this  point  should  be  placed  in  the 
middle,  that  is,  =  tn,  =  tn  +  |A t,  =  tn+1.  We  can  then  easily  write  out 

the  explicit  form  of  the  scheme 


(3.1) 


u 


(1) 

1 


u 


(2) 

1 


u 


71+1 


/  +  ^A  tL(un) 
u?  +  ^A  tL{u^) 

un  +  \a t  L(un )  +  h(u?)  -  L(u <2))) 

u?  +  T^iAt  (l(u2])  -  L(u  S1}) 
un  +  ^At  (^L(un)  +  ~L(v?)  -  ^(42))) 
41}  +  \o2At  (t(41})  -  A(41}))  +  ^A  t 


)  +  2At  +  3L(M^}  + 


For  our  analysis,  the  following  equivalent  form  of  the  scheme  is  more  convenient: 
(3.2) 

u[1]  =un  +  ^A  tL(un) 

u{2)  =  u[L)  +  -A  tL(u^) 

u^]  =un  +  ^A t  (^-L(un)  +  ^L(uS1})  -  ^l(mS2))^ 

u2  ]  —  +  ~2lA^  '))  —  L(u[  ))sj  +  -At  L(un )  +  -T(wi  +  ^(ui  ^ 

«31}  =un  +  ^At  ^L(un)  +  ^L(41})  -  Y^l(m22)) 

L(41))-L(41)))  + 

We  now  attempt  to  rewrite  each  equation  in  (3.2)  as  a  convex  combination  of 
forward  (or  backward)  Euler  steps,  as  in  the  previous  section.  The  first  two  equations 


un+l  =  un  +  l_Q2^t  ^ 


-At 


^(■U'1)  +  -T(«2  ^  +  g^(M2 
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are  already  of  the  forward  Euler  type  and  would  be  SSP  for  a  CFL  coefficient  c  =  2. 
We  would  need  to  write  the  remaining  equations  for  u^\  u£\  u A  and  un+1  into 
convex  combinations  of  forward  (or  backward)  Euler  steps.  We  present  the  details 
of  this  procedure  for  the  last  equation  involving  un+1  only,  as  the  process  is  similar 
for  the  other  equations. 

To  this  purpose  we  take 

(3.3)  >  0,  afl  >  0,  >  0,  ctg  >  0,  ctg  +  ctg  +  +  ctg  =  1, 

and  further 

(3.4)  I3‘H  >  0,  /?g  >  0,  >  0,  /f>  +  /3g  +  /jg  =  ag , 

and  rewrite  the  Erst  term  un  on  the  right  hand  side  of  the  last  equation  in  (3.2)  as 

Un  =  («3J  +  «322  +  +  «32l)«n  =  (P®  +  +  Pzl  +  «322  +  +  “S)"" 

After  a  further  algebraic  manipulation  using  all  the  equations  in  (3.2),  we  can  then 
rewrite  the  last  equation  in  (3.2)  into  the  form 


un+1  = 


+  l --  —a{2)  -  -a(2)  -  — a(2)  -  V2)  -  -B{2) 

P3,iu  +  \  g  24tt3,2  gtt3,3  24tts’4  2'°3,2  2^3,3 


(3.5) 


+ 


+ 


+ 


«322«21}  + 


(2)  (2)  . 
a3,3U2  + 


S1}  + 


Pfln  ?  + 


24 

2  1 
- 6i 

3  2 


A  tL(ur 


\9la3 


(2) 

3 


1  (2) 

3«s,4 


A  tL{u 


(W 


+ 


1 

24" 


-a 


V\a3,3  3^3,2  3^3,3  2'  3,3 

(2)  W(2A  \+t  f„,(2) 


A  tL(u 


(2) 


(2) ' 
2  , 


+ 


^(2)  ( 1)  _| _ n  A  /  T  (n  ,(■*■)  \ 

a3,4u3  +  U2^tL{U3  ) 


(2) 


?(2) 


2 

A  tL{u 


(W 

i  > 


Pag  -  )  A tL(uY') 


To  simplify  and  standardize  the  notations,  we  denote 
(3.6) 


,(i) 


..(2) 


°2,4  —  a3,2 1 


(2)  _ 
6^2  4  — 


(2) 

3,3> 


,(1) 


,,  (2) 


°3,4  —  a3,4’ 


l(0) 

'1.4 


—  3 ^ 

~  P*3,l5 


l(1) 

*1,4 


=  A 


(2) 
3,2  5 


*1,4 


=  A 


3,3 


and 


A,  4  —  2  2^2  2^ia3,3  2^3,4’ 


.,(2) 


.,(2) 


7  (2)  _  1  .  1  (2) 
P1  “  6  24“3'4’ 


(3.7) 


A1)  _  a 

°3,4  -  2^5 


6(0)  =  -  - 
1,4  6 


°1,4  ~ 


4«322  +  ^1«323 


0  (2) 
24a3,2 

2  (2) 
3«3,3  - 


_  1^,(2) _ _n^)  _  1/3W  _  l/iw 

0^3,3  24U3,4  2^3,2  2^3’3’ 


,,(2) 


}(2) 


}(2) 


-3{2) 

W  3,3  5 


&(2)  - 
°1,4  — 


1  (2) 
24a3,2 


-Q: 


(2) 

3,3 


and  write  (3.5)  as 
(3.8) 


un+1  = 


y. 

hj 


am  +bjAAtL(uj) 
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Similarly,  we  obtain 
(3.9) 


= 


y, 


(0 


i(i) 


h3 


with 

(3.10) 

where 

(3.11) 
Also 

(3.12) 

with 


,(o) 


„(i) 


°1,1  —  a2,l) 


(1) 

a\  [  =  a, 


(i) 

2,2  5 


6(0)  - 

°i,i  — 


24 


-a 


(i) 

2,2 


1  (1) 
2Q'2,3’ 


(2) 

a,  {  —  a 
hi1)  — 

ui,i  — 


(i) 

2,3) 


-a 


(i) 

2,3) 


b(2)  ~ 
°1,1 


1 

24 


«2d  >  0,  a^l  >  0,  >  0,  J  +  a^2  +  =  !• 


,  (2) 


U  o  = 


y  (4H0  +  bfyAtL(u 


(*) 


(0)  _ 

a  1,2  —  a2,l) 


(3.13)  4g  =  7  -  -  77<3  -  7T^2, 


(2) 

2,1) 

1  1 


,  ,  _  (2) 
°1,2  —  a2,2) 


,(1) 


(2)  (2) 
a  1  2  =  OL 


(2) 


„(2) 


.,(2) 


2,3) 

l(1) 


(1)  (2) 
°2,2  =  «2,4) 

1 


7AO  Q  „,(2)  A 

*fi,2  ~  g  cgl  2a2,3  3^2,4) 


(2) 


6  2 
7  (2)  _  1  ,  1  (2) 

6l’2  “  6  +  24"2’4’ 


24 


*2,4) 


/,(!)  _  n 

b2,2  ~  7^1 


4*1  >o,  422>o,  «g>o,yi>o, 


= 


where 

(3.14) 

And  finally 

(3.15) 

with 

(3.16) 

a(0)  -  3{1) 

ul,3  —  33,1 ) 

b(o)  _  j>_  _  5  (i)  _  - 

°1,3  “  0  o  O/l  “3, 2  />a3,3  o/J: 


(2)  ,  (2)  .  (2)  (2)  1 

a2,l  +  a2,2  +  a2,3  +  a2,4  ~  1- 


y  (ag^°  +  b$AtL(uf)) 


hj 


a(!)  _  73(1)  ( 2)  _  73(1) 

a  1,3  —  ^3,2)  ul,3  ~ '  ^3,3) 


24  24  ^  6 

7,(2)  _  1  W  1rv(1) 
&!,3  “  7^3,2  -  g«3,3) 


1  (1)  _  1«(1)  _  1/0(1) 
2^3,3) 


3,2 


a(1)  -  a(1) 

u2,3  —  *-*3,2 ) 
1 

l(1) 


,  (2) 


(1) 


a2,3  —  a3,3) 


Ki  =  -  3^. 


(1) 


.,(1) 


j(l) 


4,3  —  g  2^ia3,3’ 


(1) 


41  =  -- 

2, 3  24 


>  0,  >  0,  >  0, 


+  <41  +  °41  —  !) 


where 

(3.17) 
and  further 

(3.18)  /4,1  —  0)  4,2  —  0)  4,3  —  0)  4,1  +  /1, 2  +  4,3  =  a3,l- 

We  have  now  written  all  the  equations  in  (3.2)  as  convex  combinations  of  forward 
or  backward  Euler  steps,  depending  on  the  signs  of  b^\,  in  (3.8),  (3.9),  (3.12)  and 
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(3.15).  We  notice,  from  their  definitions  in  (3.7),  (3.10),  (3.13)  and  (3.16),  that 
&24,  bill ,  b^2  and  \  are  always  non-negative,  lf2\  and  b^\  are  always  non-positive, 
and  the  other  b^\  could  be  either  positive  or  negative,  at  least  a  priori.  Because 
of  our  stability  assumption  (1.5)  for  the  Euler  forward  step  (1.3)  and  the  Euler 
backward  step  (1.7),  we  would  need  to  replace  the  operator  L(u^k)  by  L{u'^k)  when 
the  corresponding  is  negative.  After  this  modification,  the  scheme  (3.2)  is  clearly 
SSP  under  the  modihed  time  step  restriction  (1.8)  with  the  choice  of  the  CFL 
coefficient 

(3.19)  c  =  max  min 

i,j,k 

subject  to  the  restrictions  (1.11),  (3.3),  (3.4),  (3.11),  (3.14),  (3.17)  and  (3.18). 

As  before,  we  optimize  the  equivalent  problem: 


(3.20) 


c  =  max  z 


subject  to  the  restrictions  (1.11),  (3.3),  (3.4),  (3.11),  (3.14),  (3.17)  and  (3.18),  and, 
for  all  the  relevant  i,  j  and  k , 


(3.21)  aft  > 

by  the  Matlab  routine  “fminicon”.  As  mentioned  before,  when  the  resulting  b^\ 


b(i)\ 

Jj,k  \ 


is  negative,  we  will  change  the  relevant  L(u^k)  by  L(u^).  The  optimal  scheme  i 
terms  of  the  CFL  coefficient  (1.8)  is  the  following 

(3.22) 

u?  =un  +  tL{un) 
u(2}  =u^  +  -A  tL(u[l]) 


,(*)■ 


m 


=  ^a[u(un  +  bi1u(AtL(un)J  + 
=  (^af^un  +  b^AtL(un)^j  + 
+  +  b[2lAtL(u[2'1) 

u g1-1  =  +  b^3AtL(un)^j  + 

+  +  b^lAtLiu^) 

un+1  =  (ai°!un  +  b^AtL(un)^j  + 

+  (a24M21]  +  b^AtL^) 


,(0)„ 


do) 


+  b\,{AtL(u[l,)j  +  (a^  +  b[z>AtL(u[z>) 

ap2wi1)  +  b^AtL^)^ 

+  +  b^AtLiu^)^ 

aP3wi1^  +  b[llAtL(u('2'1  +  b^3AtL(u[2'1) 

+  (a223M22)  +  b^AtL^) 

+  b^AtL(u([',)Sj  +  [a\ziu\z>  +  bKZ^AtL{uKZ>) 

+  (a224M22)  +  &224AtL(M22)))  +  (4>31}  +  b^AtL( 


u 


(1)' 
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with  the  coefficients  a^\  and  given  by  (3.6),  (3.7),  (3.10),  (3.13)  and  (3.16),  and 


-(*) 


(3.23) 


a 


(i) 


;;  =  o.29i2, 


a 


(i) 

2,2 


=  0.2911, 


a 


(2) 


a 


(2) 

2,3 


=  0.2453, 


f&2  =  0.2686, 
(if}  =  0.0811, 


=  0.5026, 
=  0.2457, 


=  0.1374, 
=  0.3664, 
afl  =  0.0000, 


a 


(2) 

2,2 

hi) 


=  0.0736, 


P#  =  0-1284, 


=  0.2435, 


a 


(2) 

3,3 


$1  =  0.1120,  e1  =  0.8393,  e2  =  0.7884. 


The  CFL  coefficient  for  this  scheme  is  c  =  1.2956.  Therefore,  we  have  proved  the 
following  result. 


Theorem  3.1.  The  third  order  DC  scheme  (3.22)-(3.23)  is  SSP  under  the  time  step 
restriction  (1.8)  with  the  CFL  coefficient  c  =  1.2956. 


Even  though  the  CFL  coefficient  for  the  scheme  (3.22)-(3.23)  is  reasonably  high, 
it  requires  10  evaluations  of  L  or  L.  Comparing  with  the  optimal  SSP  third  order 
Runge-Kutta  method  in  [11,  4],  which  has  a  CFL  coefficient  1  and  requires  only  3 
evaluations  of  L ,  the  third  order  SSP  DC  scheme  (3.22)-(3.23)  is  much  less  efficient. 
Of  course,  since  we  have  used  an  optimization  routine  to  obtain  the  optimal  value 
of  c,  we  cannot  guarantee  that  we  have  obtained  the  theoretical  optimal  value  of 
this  CFL  coefficient.  Theorem  3.1  provides  therefore  only  a  lower  bound  of  the  CFL 
coefficient  to  guarantee  SSP.  The  actual  DC  scheme  may  be  SSP  for  a  larger  value 
of  the  CFL  coefficient. 

If  our  objective  is  to  have  as  few  evaluations  of  L  or  L  as  possible,  we  may  require 
as  many  b^\  to  be  positive  as  possible.  A  careful  search  reveals  that  we  need  at  least 

9  evaluations  of  L  or  L  to  obtain  a  SSP  scheme.  This  leads  to  the  following  third 

order  DC  scheme 

(3.24) 

u?  =un  +  tL(un) 


uf]  =uf’  +  -A  tL(u\L>) 

,(i) 


.(i) 


,(ib 


<’  =  (4°y  +  CA^K)J  +  +  &£AiL(«HJ  +  +  bf>AtL(uf>) 

=  ^ai°2'u"  +  frig  A  tL(un)^j  +  ^a^2ni1)  +  ^lgAf^-u^)  j 

+  («122«12)  +  &S22A^(«l2)))  +  (4!  2«21}  +  4gA^(41})) 

=  (a^3un  +  b^3AtL(un)^  +  +  b^AtL^u^f'j  +  +  bf^AtL(u^) 

+  (a2lu2}  +  b^AtL^))  +  (423m22)  +  b^lAtLfup) 


do)  , 


(1) 


,(!)• 


,(2)„,(2) 


l(2) 


.(2)> 
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U 


n+1 


=  ( agu"  +  b^AtL(un) )  +  ( aiX'  +  b^AtL^) )  +  ( +  &$AiL«') 


l(0) 


(1), 


,(!)' 


,(2)  (2) 


l(2)  . 


,  (2) ' 


a<2^u^  +  frgAtZ^'ug))  +  +  ^gAil^wg))  +  +  b^AtL(u^ 

with  the  coefficients  ag  and  b^\  given  by  (3.6),  (3.7),  (3.10),  (3.13)  and  (3.16),  and 


(i)' 


(3.25) 


eng  =  0.5833, 
ctg  =  0.1650, 
42  =  0.1550, 
42  =  0.2827, 


ag  =  0.2041, 

ctg  =  0.6266, 

af\  =  0.3654, 
42  =  0.0602, 


ag  =  0.4310, 
ag  =  0.3065, 
ag  =  0.0593, 


«g  =  0.0000, 

/3g  =  0.3603, 
(2)  =  0.1652, 


a 


3,3 


d1  =  0.8990,  02  =  0.9115. 

The  CFL  coefficient  for  this  scheme  is  c  =  0.8990.  Apparently,  this  scheme  has  a 
much  smaller  CFL  coefficient  and  only  1  fewer  evaluation  of  L  or  L  than  that  of  the 
scheme  (3.22)-(3.23),  hence  is  much  less  efficient. 

As  indicated  in  the  introduction,  the  original  spectral  deferred  correction  scheme 
in  [2,  8]  corresponds  to  9\  —  02  —  1.  Within  this  subclass,  we  apply  our  optimization 
procedure  above  to  obtain  the  following  third  order  DC  scheme 

(3.26) 

41}  =un  +  ^AtL(un) 


u 


(2)  _„,(!) 


i  =u\l>  +  2AtL(Mii)) 


wg  =  ^agun  +  b^lAtL(un)j  +  +  &gAtZ/(wg)j  +  ^agwg  +  b[2\AtL(u[2'1) 

■ug  =  (4 °4n  +  b^AtL(un)J  +  +  &gALL(wg)  j 

+  («i24i2)  +  bi?2AtL(ui}))  +  (ag^  +  tgAti^ug)) 

■ug  =  ^agun  +  b^AtL(un)^  +  ^agwg  +  &gA£L(Mg))  +  ("ag^g  +  b^3AtL(u^) 
+  (ag^g  +  &gA  iL(ug))  +  (agwg  +  b(2lAtL[u{2) 


u 


n+l 


=  (aiV  +  +  (yiX’  +  b\%AtL{u\”))  +  +  6$A«L(uH 

+  (Alb1*  +  tibA'iii,11))  +  (4>f  +  4Au(«7))  +  (bib11  +  b'Aub 


r,(0) 


,(!)„.(!)  ,  ,(1) 


(A 


(1)' 


coefficients  ag. 

and  bj’k 

given  by  (3. 

6),  (3.7),  (3.10),  (3.13)  and  (3.16), 

a(11 

u2.1 

=  0.3333, 

«g 

=  0.3333, 

43 

=  0.1405, 

A 

=  0.1405, 

«s 

=  0.1977, 

4 + 

=  0.5636, 

43 

=  0.2552, 

AS 

=  0.1814, 

A 

=  0.2124, 

4+ 

=  0.1742, 

(2) 

«3,2 

=  0.1092, 

(2) 

4,3 

=  0.1961, 

=  0.0577, 

43 

=  0.0872, 

01  = 

1.0000, 

d2  =  1 

.0000. 

(3.27) 
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The  CFL  coefficient  for  this  scheme  is  c  =  1.0411.  However,  it  needs  11  evaluations 
of  L  or  L.  Therefore,  it  is  much  less  efficient  than  the  scheme  (3.22)-(3.23). 

Within  the  subclass  of  6 i  —  02  —  1,  we  can  also  explore  SSP  schemes  with  as  few 
evaluations  of  L  or  L  as  possible.  We  would  still  need  at  least  9  evaluations  of  L  or 


L  to  obtain  a  SSP  scheme,  namely  (3.24)  with  the  coefficients 
(3.6),  (3.7),  (3.10),  (3.13)  and  (3.16),  and 

afl  and  bfk  given  by 

a2\  =  0.5866, 

all  =  0-2058, 

all  =  0.4773, 

afl  =  0.0811, 

all  =  0.1170, 
(3.28)  ; 

Pfl  =  0.1298, 

all  =  0-6133> 
afl  =  0.3786, 

all  =  0-3112, 
afl  =  0.1799, 

fill  =  0.3535, 

«S  =  0.1170, 

Pf  =  0.2945, 

Pfl  =  0.0599, 

QA 

o 

o 

o 

o 

e2  =  l.oooo. 

The  CFL  coefficient  for  this  scheme  is  c  =  0.6491,  which  is  not  very  impressive. 

Finally,  we  consider  a  special  class  of  the  third  order  DC  scheme  (3.1),  in  which 
d2  =  0.  In  this  subclass,  we  do  not  need  to  evaluate  uf\  hence  this  may  lead  to  a 
scheme  with  fewer  evaluations  of  L  or  L.  After  removing  the  constraints  associated 
with  the  evaluation  of  ul  ^  and  setting  d2  =  0,  the  optimization  procedure  described 
above  yields  the  following  scheme  within  this  subclass 

(3.29) 

uf  =un  +  ^A  tL(un) 


uf  =u[1]  +  -A  tL(uf) 
uf  =  +  &i°iA  tL(unfj  + 

uf  =  (afl2un  +  bf2A  tL{unfj  + 
+  (af2uf  +  b[2}2AtL(u[2) 
un+1  =  (agun  +  bfAtL(unf)  + 
+  {afuf  +  bflAtL{uf 

with  the  coefficients  af\  and  bfk 


(apiMi^  +  &ipA£L(u^)^  +  + 

(aSui1}  +  bf2AtL(u  i1})) 

))  +  +  bf2AtL{uff} 

+  frglA  tL(uffj  +  ( affuf )  + 

))  +  («224M22)  +  b2AAtL(uf)^ 
given  by  (3.6),  (3.7),  (3.10)  and  (3.13), 


bfAtL^uf) 


b^AtL{ui 


(2b 


and 


af  =  0.3238, 
(3.30)  all  =  0.1774, 
Pf  =  0.0757, 


a(1) 

a2,2 

=  0.3237, 

afl  =  0-1264, 

(2) 
«  2,2 

=  0.2204, 

=  0.2825, 

afl  =  0.5589, 

(2) 

«3,3 

=  0.1586, 

=  0.2038, 

QA 

o 

o 

o 

o 

02  =  0.0000. 

The  CFL  coefficient  for  this  scheme  is  c  =  0.9515.  This  scheme  is  still  less  efficient 
than  the  scheme  (3.22)-(3.23),  even  though  it  has  only  8  evaluations  of  L  or  L,  2 
fewer  than  the  scheme  (3.22)-(3.23)  has. 


12 


Strong  stability  preserving  property  of  the  deferred  correction  time  discretization 


We  can  further  reduce  the  number  of  evaluations  of  L  or  L  to  7  within  this 
subclass,  yielding  the  following  scheme 

(3.31) 

=un  +  tL(un) 
uf]  =u^  +  -AtL(u^) 


=  (a[°lun  +  b^AtL(un)J  +  AfL(w^)  j  +  +  6®  A  tL(u^] 

uf  =  ^af^un  +  bflAtL(un)j  +  +  b^AtL^u^)^ 

+  (af^u^  +  b(^AtL(u[2)) )  +  (a^u?  +  b^AtL^)^ 
un+1  =  (af\un  +  b^AtL(un)^  +  +  b^AtL^^)^  +  +  bf\AtL(uf)] 

+  +  ^lAtL(41}))  +  (aS«22)  +  h2A^tL(u2])) 

with  the  coefficients  a-\.  and  l)-  \:  given  by  (3.6),  (3.7),  (3.10)  and  (3.13),  and 

=  0.5862,  =  0.2481,  ctQ  =  0.4252,  a[%  =  0.0293, 

(3.32)  =  0.1315,  a(32}  =  0.4546,  =  0.4281,  ctg  =  0.1173, 

42)  =  0.3387,  /3$  =  0.1147,  ffi  =  1.0000,  d2  =  0.0000. 

The  CFL  coefficient  for  this  scheme  is  c  =  0.7040.  This  scheme  is  slightly  less 
efficient  than  the  scheme  (3.29)-(3.30). 


4.  Fourth  order  discretization 


For  the  fourth  order  (s  =  3)  DC  time  discretization,  there  are  two  subgrid  points 
inside  the  interval  [tn,tn+1].  By  symmetry,  these  two  points  should  be  placed  at 
t W  =  tn  +  a  At  and  =  tn  +  (1  —  a)  At  respectively  for  0  <  a  <  We  will  only 

consider  the  standard  choice  of  the  Chebyshev  Gauss-Lobatto  nodes  with  a  = 
however,  see  Remark  4.2  for  the  general  case  of  arbitrary  a. 

With  the  choice  of  the  Chebyshev  Gauss-Lobatto  nodes,  we  can  easily  write  out 
the  fourth  order  DC  scheme 

U?  =un  +  ^A,i(/) 

u[2)  — u ^  +  ^A  tL(u^) 

5 

n?  =u?  +  ^AAtL(u?) 


(4.1) 
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tin  =UU  +  -At 


1  ,  fll  +  y/5 


25- 


25  —  13\/5 


-1  +  V5 


+  — ^(41})  +  60  v  l(u i2))  +  — ^l(u!3); 


«22)  =«21}  +  "If  Atffi  (L(M21])  -  ^(“l^)) 


+  -At  -^rLK)  +  -^rL(u^)  +  -  f-L(«n 


(3)  (2)  .  5  — 

tin  ’  =1i  H - 

2  2  10 


A td2  ^L(u[2))  -  L(u^)^ 


1  Aj  /-1  +  V5 


+(«")  + 


25  -  13\/5 


25  —  v7^ 


11  + 


AK;)  +  -^(<0  +  — 


M  =Mn  +  -At 


1  ,  /11  + V5 


25-  v7^ 


25  —  13\/5 


-1  +  V5 


+  — ^a(41})  +  60  v  ^(42))  +  — 


,(2)  +1)  ,  /r/  .(!)> 


= U 3  H - Af#3  I  +(«3  j)  —  L(u 


^  T  (vn\  ,  TVEr,*),  ,  Jy^_T  (v(2))  _  ^  T  (iS3h 

30  ’  30  ^  2  30  ^  2  30  ^  2 


(3)  (2)  , 

U\  =U\>  + 


5  —  a/5 


Af6>4  ^L(u^2))  -  ^(«22))) 


1,  /-1  +  V5 


L{un)  + 


25  -  13\/5 


25- v7^ 


11  +  v7^ 


L«>)  +  +  —^AAL(v>2 


uY1  =un  +  -At 


1 A ,  fn  +  Vb 


25-  v7^ 


25  —  13\/5 


-1  +  V5 


+  jV^a(41))  +  60  v  ^c42))  +  — 


■uf]  =-u^1}  +  ^-At65  (l(u 41))  -  L(m^ 


+  ±At  -^L(«»)  +  ^L(41})  +  ^(42))  -  -^(43)) 


wn+1  =«42)  +  — Y^-Atd6  (l( u f})  -  L(u^P 


+  +  ^lk«)  +  +±^i(+: 


We  can  rewrite  (4.1)  into  an  equivalent  form  similar  to  (3.2),  then  attempt  to 
rewrite  each  equation  as  a  convex  combination  of  forward  (or  backward)  Euler  steps, 
as  in  the  previous  section.  The  first  three  equations  are  already  of  the  forward  Euler 
type  and  would  be  SSP  for  a  CFL  coefficient  c  >  2.  As  to  the  fourth  equation,  we 
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can  rewrite  it  as 
(4.2) 


411  =  E  (4:i“f + b^tL(Uf)) 


h3 


with 


,(o) 


„(i) 


°1,1  ~  a2,l> 


(1) 

a\  [  =  a 


(i) 

2,2; 


(4.3) 


,  (o)  _  11  +  _  5 _ 

1,1  120  10 

(1)  _  25  -  y/5  VE 
Vi  — 


VE  (i)  5 

«2,2  -  - 


(2)  =  «(1) 

V5 


°1,1  ~  “2,3; 


,(3) 


.,(1) 


°1,1  —  a2,4) 


10 


- 


5  —  VE 

"To" 


b{2)  ~ 
°1,1  — 


120 

25  -  13VE 
120 


5 


v/5  (1) 

—“2,4; 


5  —  VE  (i) 

Of  9  A  « 

10  2’4’ 


b{3)  - 

°i,i  — 


1  +  V5 
~120 


where 


(4.4)  (X2J  —  0;  a2,2  —  0;  a2,3  —  0;  “2,4  —  0; 


t(1) 

*2,1 


Olo  1  +  CX2  2  T  OL^  3  +  0(94  —  1. 


1(1) 

*2,3 


*2,4 


Similarly,  the  fifth  equation  in  (4.1)  can  be  rewritten  as 


(4.5) 

with 


(2) 

u)2  = 


Y  [a%uf + h{jl^tL(u. 


(o 


/<> 


h3 


,  ,  _  (2) 
®  1,2  ~  a2,l; 


l(0) 

'1,2 

1,(0) 


,(!) 


„(2) 


°1,2  —  a2,2i 


,  (2) 


11-^5  5-VE  (2) 
“2,2 


(4.6) 


/hUj  = 

1,2  120 
(!)  _  25  +  13y/5  VE 
5  * 


CL- ^  2  — 

5- V5 


(2) 

2,3; 


10 


-« 


(2) 

2,3 


(3)  (2) 

°1,2  =  «2,4; 

5  ~  VE  (2)  _ 
10  2’4 


r« 

s2,2 

H  +  a/5 


a2,2  —  a. 


(2) 

2,5 


(2) 

120  a2’5’ 


V5 


^1,2  —  777  “T^l  1 


b{2)  ~ 

ul,2  ~ 

bi3)  ~ 

u  1,2  — 


120 

25  +  a/5  5  —  VE 

120 

i  +  VE 


(2) 

10  “2’4 

-1  +  V5 


(2) 

5  ~2’3 

25  -  13\/5 


VE  (2)  25  VE  (2) 

— “2,4  120  a2’5’ 


120 


-ct 


(2) 

2,5; 


120 


120 


-a 


(2) 

2,5; 


JO  _ 

0-2,2  ~  ~0\ 


where 

(4.7) 


V2}  >  0,  ctg  >  0,  ag  >  0,  ctg  >  0,  >  0, 


t(2) 

*2,1 


G4  i  +  Ct^  T  ®)  .‘!  +  fl)  i  +  Qi)  (  —  1. 


t(2) 

*2,3 


t(2) 

*2,4 


t(2) 

*2,5 


The  sixth  equation  in  (4.1)  can  be  rewritten  as 
(4.8)  u{3)  =  Y  {a%uf  +  b%^tLV{j 


1,3 
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with 

(4.9) 

al03  =  “S, 

(1)  (3) 

Oi,3  —  a2,2i 

(2)  (3) 

«i, 3  =  «2,3. 

(3)  (3) 

°i,3  ~  02  4, 

(1)  (3) 

a2,3  —  02)5, 

(2) 

02,3 

/,<»2  =  1  - 

5-  VE  (3) 

- ryX  L  — 

5  -  ^  (3) 

5 -VE  (3) 

- ()  \  i  — 

11  +  a/5  (3) 

- rvA  l 

11 

'  V5  (3) 
- «2,6> 


h ^ 
"1,3 

&(2) 

"1,3 

6(3) 

"1,3 


5 

12 

5 

12 

1 

12 


^  _  hLm 


25 


5 

5  —  VE 


^4%  +  #9l0w 


-0. 


(3) 

5-«2,3  g  «2,4  12Q  -*,o  ■  g  — *,o 

5  —  y/E  (3)  25  —  13\/5  (3)  25  +  -\/5  (3) 

- a n  4 - c - Clin  r. 

2,4  120  2,5  12n  2,6 


12U 

25  +  13^  (3) 
- — - <x*  L 


1L20  “2’6’ 


where 

(4.10) 

u2,l  T  «2. 
The  seventh  equation 

(4.11) 


10  "  10  2)4  120  2" 

~  1  +  Aj3)  ,  1  +  ^n.O)  /  (1)  _  [n  _  A)  (3)  ,  (2)  _  5~  \ 

120  2’5  +  120  2’6’  &2’3  “  5  0l  5  ^  2’6’  &2’3  “  10 


v/5 


-0, 


5  >  0,  ag  >  >  0,  >  0,  cxg  >  0,  ag  >  0, 

43+^  =  i- 


«2J  >  0,  <*2,2  >  0,  «2,3  >  0,  C 

(3)  ,  (3)  .  (3)  .  (3)  .  (3)  .  (3) 

a2,l  +  a2,2  +  0-2,3  +  02)4  +  a2,5  +  a2,6  : 

equation  in  (4.1)  can  be  rewritten  as 

12  (aJAuJZ)  +  b%ML{uf 


A 


with 

(4.12) 

(0) 

<4  = 


(0)  _  o(l) 

M3,l  5 


(2)  (1) 
4,4  =  <*3,3, 


Jl)  _  M  1)  _(2)  _  M 1)  „(3)  _  Ml)  (1)  _  (1) 

/"3,1,  ul,4  —  M3, 2)  “i,4  M3, 3|  1*1,4  —  M3, 4,  u'2,4  —  <*3,2  > 

(0)  11  +  VE  11  +  y/E  (1)  11  —  y/E 

120  120  “3’2  120  03,3  12~d’4  10  ^  10 

25  +  13^  (1)  ,  VEn  (d  5  a)  ^ 
120  Q'3’3  +  5  ia3’4  12%’4  5 

c  .  /E  c.  c  .  fE 


"1,4 


j/1) 

"1,4 


120 

25  -TLm  +  //agj 


1  (i) 

rf3A 


(3)  (1) 

4,4  =  <*3,4 

5  -  ^  «  5  -  VE  (i) 

M3, 3  M3, 4  5 


b(2) 

"1,4 

b(3) 

"1,4 

"2,4 


25 


120  “3’2  '  5  i: 

13 1/5  (i)  25  +  y/E  (i) 

- CX o  9  —  - CXo  q 

20  3’2  120  3’3 


10 

5  ,(I)  _  VSjg  _ 


3,4  , 


10 


120 

1  +  a/5  (!)  1  +  a/5  (!)  1  (!) 

<*3,2  +  a3,3  jj®3,4’ 

25  -  13\/5 


>  _.(i)  5  _.(i)  5  +%(*) 

-"2<*3,4  ^OL 3i4  ^  M3, 4 , 


120 

25-  V5 


120 


5 


+i<4':5  ~  ^i<*£L 


6(2) 

"2,4 


5  ^41 


A3) 

"2,4 


where 

(4.13)  ag  >  0,  >  0, 


d1)  >  n  >  0 

-3,3  —  u>  <*3,4  —  u> 

and  further 

(4.14)  >  0,  /3$  >  0,  /3$  >  0,  PI]A  >  0, 


120  10 

°4q  +  0^2  +  0^3  +  a3A  =  1) 

0$  +  ^2  +  43  +  0$ 


■ l  +  \ / 
120 


a{1) 

a3,l- 
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The  eighth  equation 
(4.15) 


(2) 

«3 


■Ei 

hi 


+  b^AtL(u{-- 


with 

(4.16) 

(o) 

°1,5 

a(1) 

u2,5 

6(°) 

°1,5 


b (1) 

°1,5 


6(2) 

°1,5 

6(3) 

°1,5 

6(1) 

°2,5 

6(2) 

°2,5 


_P3,l5 


(2) 

=“3,25 

11 


a(1)  -  /?(2) 

ul,5  —  ^3,2  5 

7(2)  _  (2) 

'  —  “3,3  5 


(2) 

«!, 5  = 
(3) 

a2,5  ‘ 

11 


25  a2,5  ^3,35 

—  VE  n  +  VE  (2) 

120  120  “3’2 

2) 


-  /?(2) 

‘  M3, 3  5 

(2) 

=  «3,45 

'  V5  (2) 

- “3,3 


(3) 

«1,5 


■  /3(2) 

‘  /A, 4  5 


,M  -  a(2) 

"  -  ~  “3,5 


5 


'  V^fl(2) 
10  ^3’2 
25-^5  (2) 


10 

A, 


120 

5  _  v^5„(2) 

'M3, 3  -^Q  M3, 4  5 


«3,5 

1  (2) 
-  ^3,4 


11  +  VE  (2) 
120  “3’5 


,,(2) 


120^  +  5oA,3 

5  ,(2)  _  ^fl(2)  _  y^_p{ 2) 

■3,4  g  “3,3  5  “3,4 


10  ' 

25  +  W5  (2)  ,  VEn  (2) 
120  3’3  5  3’4 


12«3,4 


25 


,4  g  >6,6  g  5 

13\/5  m)  25  +  a/5  m)  5  —  VE  /2 

120  3’2  120  3’3  10  3’ 
^  ^  1 

£^,(2) 

12 


120  ^  120 
■  i  +  VE  (2)  i  +  VE  (2)  w 
120  3’2  +  120  a;3’3  ^*3,4, 


A.  -,(2)  _  5  42) 

",4  12  3’4 


5  —  \/5 

10 


120  rzu 

25  +  13\/5  V5n  V5„  (2)  (2) 

"  120  5^3  "  T^‘“3’3  "  ^‘“3’4  ” 

25  +  VE  5-y/5  (2)  25  -  13\/5  (2) 

= - Qo a\  i - ol,  L 

120  10  3’ 


25 


120 


VE  (2) 

— “3,55 


b(3)  -  - 

°2,5  — 

where 
(4.17) 
and  further 


)  \  0  O  \  O  ^2)  lo\/0  /2) 

"120  lo- 1 92“3’4  120  a3’5’ 

1  +  Vjj  _  ~  1  +  vj-,(2)  7 (i)  =  A 

120  120  3i5’  3>5  5  °3 


“S  >  0,  «g  >  0,  ctg  >  0,  ctg  >  0,  afl  >  0, 

“S  +  “S  +  “323  +  “S  +  “S  =  X5 


cuivi  1  ur  uiroi 

(4.18)  $1  >  0,  422}  >  0,  / i2l  >  0,  42i  >  0, 

The  ninth  equation  in  (4.1)  can  be  rewritten  a 

(4.19)  uf*  =  V  (afyuf  +  b{-l&tL(uf 


Y.  liu,  C.-W.  Shu  and  M.  Zhang 


with 


(4.20) 

a(0) 
a  1,6 

~d{3) 

— +  3,ii 

a(1)  -  d(3) 

“1,6  —  +3,2  > 

a(2)  -  3{3) 

“1,6  —  +3,3  > 

(3) 

«1,6 

—  /3(3^ 

“  +3,4  > 

a(1) 

u2,6 

(3) 

=«3,2» 

a(2)  -  a(3) 

u2,6  —  “3,3  > 

45  =  45  > 

a(1) 

u3,6 

(3) 

=  <*3,5  > 

ul,6 

1 

11  +  ^  (3) 

ii  -  VE  (3)  l 

n(3)  - 

11  +  v7^ 

“12  ~~ 

120  3’2 

120  3’3  12 

“3,4 

120 

~a3,5 


«(2)  -  a(3) 

u3,6  —  “3,6 
11  -  V7^ 


5  ^+3)  5  ^0(3)  5  ^o(3) 

A-^3,2  n  Y3,3  n  +3,4 


120 


-a 


h<4) 

°1,6 


6(2) 

°1,6 

h( 3) 

°1,6 


10 

25  -  V5  (3) 


10 

a/5 

v  /I 


120 

^+3)  _  ^+3) 
5  +3,3  p;  +3,4 


(3) 


~a3,2  +  ~^la333 


10  ' 

25  + 13^  (3)  VE  (3)  5  (3) 

120  “3’3  +  5  6ia3A  12“3’4 


25 


>,0  g  /  0,<i? 

—  13V5  (3)  25  +  y/E  (3)  5  —  a/5  (3) 

T20  “3’2  120  “3’3  10  2“3’4 

1  +  VE  (3)  1  +  a/5  (3)  1  (3) 

"T7TZ - °+?  H - -  ™  a3,3  _  ^2^3, 4’ 


5  (3) 

12Q'3’4 


5  -  V 
10 


120 


6(1)  — _ 

&2'6“12 


120 

25  +  13^  (3) 

- — - «3,6, 


,(3)  _ 
",4 


25 


VE  (3) 

- «3,5 


b(2)  ~— 
°2’6“12 


/5  (3) 

— 03«3,6  - 

5  -  y/E 


-Oa 


120 
5  —  VE 


b(3)  — 
°2’6“12 


b(2) 

°3,6 


5  —  VE 


120 


10  10 

1  —  1  +  VE  (3)  1  +  VE  (3) 

a3,5  +  12O  a3,6> 


^2«33] 


25 


120 


13^  (S)  25  +  VE  (3) 

- a\  k - CKo  L 

20  3’5  120  3,6 


120 


}SX)  _  v  "a 

°3,6  -  —U3 


±ZjU 


where 

(4.21) 


10 


4, 


$  >  0.  ag  >  0 


ag  >  0,  ag  >  0.  ag  >  0.  ag  >  0.  ag 
ag  +  ag  +  ag  +  ag  +  ag  +  ag  =  1. 

and  further 

(4.22)  45  >  0,  45  >  0,  45  >  0,  45  >  o,  45 + 45 + 45 

The  tenth  equation  in  (4.1)  can  be  rewritten  as 

(4.23)  u4]  =  {a^juf1  +  b^j/S.tL(u^ 
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with 

(4.24) 

a(0)  - 

ttl,7  — 


■V(1) 

■74.11 

q  (1)  -  -d) 
ul,7  —  74,2i 

o(2) 

ttl,7 

=  773i 

(3)  (1) 

«1,7  =  74,4, 

a(1)  - 

u2,7  — 

■  R (1) 

74,3i 

a(3)  -  d(1) 

u2,7  —  74,4 1 

a(1) 

u3,7 

=  til 

n(2)  -  «(1) 
a3,7  —  *0-4,3, 

n{3)  - 
tt3,7  — 

11  + v7^ 

11  +  V5  (D 

«4.2 

11 

-V5  (i) 

.  ™  <*4.3 

V1)  11  +  v/V) 

1  ~  0*4.4  1  ™  74.2 

>n 


1  /?(!)  5  ~  (i)  5  —  y/E  (i)  5  —  v7^  (i) 

12^4,4  10  74,2  10  74,3  10  74,4  ’ 

25  —  ci')  \/5  cii  25  4“  13v^5  / 1 \  \/5  cd  5  cd 

- —Pa!  +  —OiPil - —Pi!  +  —OiAl - til 

120  4’2  5  4,3  120  4,3  5  4,4  12  4,4 

_  v7^  (i)  _  y/5  (i) 

T07’3  T^4’4’ 

25  -  13V5  (1)  25  +  V5  (1)  5  -  y/E  (1)  5  (1)  5  -  y/E__(1) 

120  74’2  120  74’3+  10  1 2  ‘4’4  To  74,45 


6(3)  -  - 
°2,7  “ 


VO  (1)  _  vo  (1) 

T07’3  V74’4’ 

25  —  13\/5  cii  25  +  y/Z  cii  5-^  m  5 

-  120  ^  +  —i r ^  ~  t/‘a 

-  1  +  (!)  1  +  y/5  (!)  1  (!) 

120  74’2+  120  74’3  12/4’4’ 

25  —  \/5  (1)  ,  y/En  (i)  25  +  13\/5  (i)  ,  y/En  (1) 

120  Q'4’2  +  V 3" 4’3  120  ',|3  +  V^4’4 

25-1325  (1)  25  +  25  (1)  5  -  25  (1)  5  (1) 

120  4’2  “  120  13  10'““  “  I?4'1 

-  1  +  a/5  (!)  1  +  a/5  (!)  1  (!) 

120  4’2  120  4’3  12  4’4’ 

-  V5  ti  (!)  y/5  (i) 

i20 —  V"3“«  "  T^“4'4’ 


72)  _25 

-  13y/5 

5  -  V5n  _  (i) 

2()  $40+41 

6(3)  -  — 

1  + 

°3,7  “ 

120 

°3,7  — 

120 

where 

(4.25) 

til  >  0, 

0*4,2  —  0,  <*4,3  >  0, 

o' 

Al 

3 

O4J  +  Q-il  d"  O4J  +  CX-il  =  1, 

and 

(4.26)  ft}  >  0,  p[]l  >  0,  /3$  >  0,  p[]l  >  0, 

and  further 

(4.27)  ji]l  >  0,  77)  >  0,  7g  >  0,  7?i  >  0, 


P$  +  P$+P$  +  p£i  =  aW, 
747  +  7^2  +  til  +  7^4  =  til 
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The  eleventh  equation  in  (4.1)  can  be  rewritten  as 
(4.28)  42)  =  51  (4X4)  +  bifiAtL(uj‘ 


with 

(4.29) 

(0)  (2) 
<8  =74,1* 

,(!)  _o(2) 
'2,8  —  Y4,2) 


fl2,8 

a(2) 

u3,8 


(1)  (2) 
<8  =  74,2, 

.(2)  _  M) 

'  Y4,3  5 


(2)  -a(2) 
■a4,3; 


a2,8 

a(3) 

u3,8 


(2) 

^4,4, 


6(°) 

°1,8 


A1) 

°1,8 


11 


il  “3, 8  4,4  5 

—  y/E  11  +  V^5  (2) 

120  “4-2 

n-A,P> 


ij 


(2)  (2) 
<*1,8  =  74,3 , 

a(3)  -  b(2) 

u2,8  —  Y4,4> 

(1)  _  (2) 
a4,8  —  Ct4,5 


(3)  (2) 

<*1,8  =  74,4, 

(1)  _  (2) 
a3,8  ~  a4,2; 


120 


_L.  W  -L  Z-lW  JL 

H  +  V^  (2)  11 -a/5  (2) 

120  /4’2  120  /4’3 


11  —  y/E  (2)  1  (2)  11  +  a/5  (2) 

7,2  120  "4’3  ~~  12“4’4  120  "4’5 

5  —  y/E  (2)  5  —  y/E  (2)  5  —  y/E  /2) 

_ ^/v  '  _  _ r\ix  '  _  _ r\i v  ' 

1 4,3  j^q  / 4,4  ? 


— /$  + 


>  _  _  I _ _  : 

!  1274,4  10  74,2 

25  +  13V5  (2)  a/5  «(2) 

120  ^4,3  +  -jMv 


a/5  (2)  a/5  (2) 


10 


6(2)  -  - 
°1,8  “ 

6(3)  -  - 

fh,8  — 

j/1)  . 

°2,8  ' 


4  g  /  4,3  g  /  4,4  * 

!3 y^5  (2)  _  25  +  a/5  (2)  5  -  yo  (2)  _  5  (2)  _  5  -  a/5j2) 

20  74,2  120  4,3  +  10  274,4  1274,4 


10 


120  '  120  ^  10 

1  +  a/5  (2)  1  +  a/5  (2)  1  (2) 

120  /  4,2  +  120  74,3  12/4,4’ 

25 -a/5  (2)  a/5  (2)  25  +  13^5  (2)  a/5  (2)  5  (2) 

“4,2  +  “30:4,3 - ^ - ft  1,3  +  “30:4,4  - 


2 

-74,4* 


120 

_  PL  A2)  _  PL  «(2) 

g  ^1P4,3  g  P4, 4  5 

(2)  25  -  13\/5  (2)  25  +  V5  (2)  5  -  y/E  (2) 

°2,8  —  220  °4’2  120  Q"4’3  '  in  y4a4,4 

(3)  —  1  +  y/E  to\  1  +  y/E  to\  1 


°2,8 

°3,8 

b(2) 

w3,8 

6(3) 

°3,8 


where 

(4.30) 


(2) 
h44,4 , 


120  4,2  120  4 

“  1  +  V^..(2)  ,  1  +  Vd  (2) 

120  4,2  +  120  4,3  12"  ;-  ‘- 

25  +  1 3-^/5  y/En  y/En  (2)  (5 

=  120  5~*5  ~  ~03"4,3  "  V^4’ 

r,  r'  1  /P"  r-  fr  r\  1—  -10  /P" 


5  (2)  _  5  —  y/E  (2) 

22  4,4  2  ( )  r72t44i4 , 


120  5  “  5 

25  +  ^  5-v^  (2)  25 

120  10  ^ 94"4,4  ~~ 

1  +  \/5  —  1  +  y/E  (2) 

- a/\  L. 


13a/5  (2) 

- «4,5> 


(2)  25  v/5  (2) 

'  :,4  120  a4,5’ 


120 


120 


a4,5., 


120 

h(1) 

w4,8 


A 


"$5 


ag  >  0,  >  0,  o%l  >  0,  af}  >  0,  ag  >  0, 

(2)  (2)  (2)  (2)  (2)  1 

&  4,1  +  Q^4}2  +  ^4,3  +  ^4,4  +  ^4,5  —  1? 
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and 

(4.31)  pfl  >  0,  pfl  >  0,  pfl  >  0,  42> 

and  further 


(4.32)  7g  >  0,  7g  >  0,  7g  >  0,  >  0,  7g  + 

Finally  the  twelfth  equation  in  (4.1)  can  be  rewritten  as 

(4.33)  un+l  =  iKa%Vj{f>  +  bj$AtL(uj 


(2)  .  (2)  .  (2)  (2)  _  fl(2) 
7,1  +  74,2  +  74,3  +  74,4  —  04,1  • 


with 


(0) 

1,9 

(3) 

=74,1, 

a(1) 

ul,9 

(3) 

=  74,2, 

(2) 

<9 

(3) 

=  74,3, 

(3) 

<4,9 

(3) 

=  74,4, 

a(1)  ■ 

a2,9  ‘ 

-  d(3) 

—  74,2 , 

(2) 

a2,9  — 

(3) 

2,9 

— gW 

— M4,4  5 

a(1) 

u3,9 

(3) 

^4,2’ 

(2) 

«3,9 

(3) 

=  «4,3, 

(3) 

°3,9 

(3) 

^4,4’ 

a(1) 

u4,9 

(3) 

=  «4,5, 

(2) 

°4,9  = 

6(0)  -  — 
6l-9_12 


11  —  V5  (3) 

CXa  o  OLa  q  A  ^  a 

120  4’2  120  4’3  12  44  120  44  120  4)6 

nW5  (3)  _  11  -  a/5  (3)  _  1  (3)  _  5  -  n/5  (3)  _  5  -  vg  (3)  _  5  -  \/5  (3) 

120  4’2  120  4’3  12^4’4  10  74’2  10  74’3  10  744  ’ 

^2^45  +  47,45  - 25  +  13v^45  +  77,45  _  I45 

120  4,2  5  4,3  120  4,3  5  4,4  12  4,4 


1  (3) 

— rv\  .  — 


11  -V5  (3) 

120  "4’6 


fe(3)  _  _ 
ui  ,9  — 


:3)  _  7LJ3) 

1,4  ]^2a4,4 


6(3)  -  - 
°2,9  — 

,(1)  _  5 

^3,9-12 

b(2)  ~  — 
°3'9  ~  12 


7  (3)  _  1 
63,9-12 


_  _ _ I  v  7/)  M3)  _  —  '  —  v-7(3)  V  _/)  M3)  _  _7(3) 

120  74’2  +  5  '  14  120  74’3  +  5  1  4-4  12  ;  ‘ 

__  V5  (3)  _  V7^  (3) 

T^4’3  T^4’4’ 

25  ~  13 \/5  (3)  _  25  +  a/5  (3)  5  -  ydj  (3)  _  J>M3)  _  5  ~  (3) 

120  74,2  120  4,3  +  10  2744  1274’4  10  74,4’ 

_1  +  V/5  (3)  i  +  V7^  (3)  1  .(3) 

120  /  4’2+  120  /  4’3  12  ; 

25  —  V5  (3)  VS  (3)  25  + 13-\/5  (3)  \/5  (3)  5  (3) 

- TbcT^'4’2  +  ip3"43 - 120 — “43  +  IT03"4-4  “  12"44 

_  PL  fl(3)  _  PL  73(3) 

g  “lP4,3  g  ^lP4,4  5 

25  -  13a/5  (3)  25  +  n/5  (3)  5  -  VS  (3)  5  (3)  5  -  VS  (3) 

120  4’2  ~  120  4’3  +  ~  T^4  -  —^-0204,4, 

-  1  +  v7^  (3)  1  +  a/5  (3)  1  (3) 

- cb  9  H - cti  4 - a  a  At 

120  4’2  120  4’3  12  4’4’ 

5  \/5  (3)  \/5  (3)  25  -  v/5  73)  Vp.  (3)  25  + 13\/5  (3) 

12  5  05  5  03  4’3  5  03  44  120  44  +  5  05  4’6  120  4’6’ 

5  5  —  Vs  5  —  Vs  (3)  25-13^  (3)  25  +  ^  (3) 

- c/6 - 67  ad  f - ad  t - ad  L 

12  10  10  44  120  44  120  4’6’ 

1  -  1  +  V7^  (3)  1  +  v7^  (3)  (l)  __  VS  VS  (3)  ,  (2)  5  -  Vs 

12  120- Q44  +  ^O-044’  64,9  -  —  05  -  —  ^4,6,  64.9  —  IQ 


5  -  A,  o(3) 

2Q 
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where 

(4.35) 

and 


a 2  >  0,  ag  >  °>  ^  °>  >  °>  «2  >  0,  «2  >  0, 


„(3) 

*4,1 


QLa  i  +  a4,2  +  0(4  s  +  CI4  ,(  +  Ot\  r\  7~  Cx\  a  —  1, 


73) 

*4,3 


73) 

*4,4 


1(3) 

*4,5 


73)  - 

*4,6 


(4.36)  (3^1  >  0,  (3^2  >  0,  >  0,  /3$  >  0,  (3^J  +  +  Ai^s  +  A47  —  ex 

and  further 

(4.37)  7®  >  0,  772  >  0,  7^  >  0,  7^j  >  0,  7 44  +  7$  +  7$  +  7 2  =  Alp- 

Similar  to  the  third  order  case,  we  can  formulate  the  optimization  problem  (3.20), 
subject  to  the  restriction  (1.11),  (4.4),  (4.7),  (4.10),  (4.13),  (4.14),  (4.17),  (4.18), 

(4.21) ,  (4.22),  (4.25),  (4.26),  (4.27),  (4.30),  (4.31),  (4.32),  (4.35),  (4.36),  (4.37),  and 

(3.21) ,  and  solve  it  using  the  Matlab  routine  “fminicon”.  We  obtain  the  following 
optimal  scheme 

(4.38) 

uP  =un  +  5  ~y^ A tL(un) 


?(3) 


?(3) 


?(3) 


?(3) 


?(3) 


j(3) 


?(3) 


.,(3) 


'•  1 


10 

s/5. 


<42)  =um  +  UiA  tL(u?) 
5 

.„(3)  =um  + 


u 


1 

(1) 


10 


^  =  (af}un  +  b^>AtL(un))  +  +  6$  AtL(u[L>))  +  +  b^AtL^ 


7°) . 


(i) 


,(i) 


72)„, (2) 


l(2)  , 


,(2) 


+  (a^  +  &i>fL(^)) 


u{?  =  (o^un  +  b^>AtL(un))  +  +  b\%AtL{u?>))  +  ( a\^>  +  b\^AtL{u\ 


7°)„ 


do) , 


(1) 


.(1) 


,(2)f2) 


l(2) 


.(2) 


+  (  al32Ml3)  +  )  +  (a2,2U2’  +  b{2^AtL(uK2 


,(3) 


,aui)  /.(l)  . 


,(!) 


U 


f  =  la^>un  +  b^>AtL(un))  +  (aiX;  +  6SAiL«0)  +  (  +  ^AfL(^ 


7°). 


7°) . 


71)„71)  ,  A1) 


(i) 


,(2).. (2)  .  j (2) 


d2) 


+  (ai33Mi3)  +  bi33AtL(ui3)))  +  +  b2^tL(u2}))  +  (a223M22)  +  b^AtL^ 


(2) 


u 


?»  =  Ui>”  +  6$AiL(u"))  +  (aS>r  +  5i;iAiL(«H)  +  («',>«+  b^tHn'i 


70). 


70)  . 


71)»71) 


(i) 


7b 


72).72) 


72), 


72) 


+  (ai34ui3)  +  b^AtL(u  (i’))  +  [a2^u2L>  +  b2’AAtL{uK2>))  +  [a2tiu2’  +  b2^AtL(u2 


1,(3) 


73) 


71)„71) 


(i). 


7b 


72),  72) 


72) . 


72) 


Q234M23)  +  ^27  AtL(l43))) 
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U3]  =  (4?5W"  +  bii°lAtL(un))  + 
+  (<$«?’  +  fcgA  tL(u? 

+  («235m23)  +  b^AtL^ 


u3  —  [a\,6un  +  b\fiAtL(un)J  + 

+  (agwg  +  6gAiL(ug 
+  (aifiu23)  +  bi3eAtL(u(23) 
u g  =  (a[jUn  +  b^AtL(un)^j  + 
+  fagwg  +  &gAtL(ng 


+  ^ a^u g  +  b^jAtL(u^ 

+  +  bi'3jAtL(uf'> 

( „(°)».«  _l.  h(°)  A +7 i 


«SMS1}  +  ^gATLl^i1)))  + 
I  +  (ag^g  +  b^lAtLiu^ 
I  +  ( as, Ws]  +  b^AtL^ 
agwg  +  6gA£L(ug))  + 
I  +  (aggg  +  b^Atl^u^ 
1  +  (aggg  +  b^AtLiu^ 
aijug  +  &gA77(-ug)  j  + 
I  +  ^agug  +  b^AtLiu^p 
I  +  ^agtig  +  b^jAtL(u^ 


a%u?  +  b%AtL(u?)) 

I  +  (ag-ug  +  b^5AtL(u{2 


a\:>uY>  +  b\:>AtL(uY>)) 

I  +  ( agw g  +  b^6AtL(u(22) 
1  +  (a3fiu3]  +  bflAtL{uf] 
a'ljU'2'1  +  &gAtL(«g)  j 
I  +  ^agug  +  frgAtl/^g 
I  +  (afjU^  +  bfkiAtL(uf'> 


=  (ai>n  +  ^AtL(nn)J  + 
+  ^ag-ug  +  6®  AiL(ug 
+  {a%su2)  +  b{2lAtL(u{23) 
+  («!343)  +  b^8AtL(u{3) 
un+l  =  (ag«n  +  &gAtZ(un))  + 

+  (4g«g  +  &i3gAtL(fff) 
+  (°Sm23)  +  b{3lAtL{u{3) 
+  (aSM33)  +  bflAtLiuf 

with  the  coefficients  ag  and  &g 
(4.24),  (4.29),  (4.34),  and 


^>^  +  ^AtL«)J  + 
)  +  (oggg  +  b^AtL(u{2] 
)  +  (aggg  +  b^AtL(u^ 
)  +  (oggg  +  b^8AtL(u^ 
(aggg  +&SAtZ(«!1)))  + 
)  +  (a£gg  +  b{2lAtL(u{2] 
)  +  (4>31}  +  4!9AtZ(41} 
)  +  (“ggg  +  &gA£L(ug 

given  by  (4.3),  (4.6),  (4.9), 


a^u?>  +  b\%AtL(uY>)) 

I  +  ( agu g  +  frgAtZ^g 
I  +  (4?8m  g  +  bflAtL(u g 


{a\:>u?>  +  b\:>AtL(uY>)) 

)  +  (aggg  +  42<WZ(wg 
)  +  (aggg  +  bflAtLiuf 
)  +  (a?9“i2)  +  b^gAtL(u(2) 
(4.12),  (4.16),  (4.20), 


(4.39) 


«2g 

=  0.2505, 

«g  =  0.2506, 

=  0.2505, 

(2) 
a  2,i 

=  0.1288, 

«g2 

=  0.1206, 

ctg  =  0.2805, 

«g 

=  0.0650, 

(3) 

«2,1 

=  0.0721, 

«2g 

=  0.1273, 

ag  =  0.0642, 

ag 

=  0.1162, 

=  0.2689, 

=  0.5507, 

ag  =  0.1061, 

«g 

=  0.2015, 

/3g 

=  0.0751, 

=  0.2270, 

gg  =  0.1554, 

=  0.2758, 

<43 

=  0.0243, 

«S 

=  0.2882, 

ag  =  0.0388, 

=  0.0666, 

43 

=  0.1137, 
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=  0.0901, 

(3) 
a  3,1 

=  0.1384, 

=  0.0464, 

43 

=  0.0160, 

“m 

=  0.1188, 

(3) 

«3,5 

=  0.2048, 

43 

=  0.0312, 

43 

=  0.0589, 

43 

=  0.0356, 

a(1) 

a4,l 

=  0.7232, 

a(1) 

a4,2 

=  0.2164, 

«(1) 

a4,3 

=  0.0473, 

43 

=  0.1778, 

43 

=  0.4305, 

R (P 

P  4,3 

=  0.0348, 

'yW 

74.1 

=  0.0224, 

43 

=  0.1321, 

43 

=  0.0104, 

43 

=  0.3139, 

(2) 

«4,2 

=  0.0000, 

(2) 
a  4,3 

=  0.2877, 

(2) 

aA,A 

=  0.0386, 

43 

=  0.1169, 

Pa, 2 

=  0.1002, 

d(2) 
Pa, 3 

=  0.0910, 

74,1 

=  0.0387, 

(2) 

74,2 

=  0.0539, 

(2) 

74,3 

=  0.0232, 

(3) 
«  4,1 

=  0.2040, 

43 

=  0.0315, 

(3) 
«  4,3 

=  0.0673, 

(3) 

aA,A 

=  0.1130, 

(3) 

«4,5 

=  0.2507, 

43 

=  0.0691, 

Pa, 2 

=  0.0653, 

g (3) 

74,3 

=  0.0595, 

(3) 

74,1 

=  0.0168, 

43 

=  0.0362, 

(3) 

74,3 

=  0.0160, 

6l  = 

0.7043, 

9-2  = 

:  1.0000, 

93  =  0.66  22, 

d4  =  1.0000,  e5 

=  0.6388,  96 

The  CFL  coefficient  for  this  scheme  is  c  =  1.2592.  Therefore,  we  have  proved  the 
following  result. 


Theorem  4.1.  The  fourth  order  DC  scheme  (4-38)-(4-39)  is  SSP  under  the  time 
step  restriction  (1-8)  with  the  CFL  coefficient  c  =  1.2592. 


The  optimal  scheme  (4.38)-(4.39)  needs  21  evaluations  of  L  or  L.  We  can  also 
obtain  a  fourth  order  SSP  DC  scheme  with  19  evaluations  of  L  or  L ,  however  the  CFL 
coefficient  is  only  c  =  0.6775,  hence  it  is  much  less  efficient  than  the  scheme  (4.38)- 
(4.39).  For  the  original  spectral  deferred  correction  scheme  in  [2,  8]  corresponding  to 
9i  =  9-2  =  63  =  64  =  65  —  9q  —  1,  we  can  obtain  a  SSP  scheme  (4.38)  with  different 
choices  of  parameters  than  those  in  (4.39),  with  21  evaluations  of  L  or  L  and  a 
CFL  coefficient  of  c  =  0.9463.  This  is  again  much  less  efficient  than  the  scheme 
(4.38)-(4.39).  We  do  not  list  the  details  of  these  schemes  to  save  space.  Finally, 
when  95  =  96  =  0,  we  do  not  need  to  evaluate  u 4^  and  u%\  leading  to  the  following 
scheme 


(4.40) 

(i) 


u 


=un  + 


u 


(2)  -A1) 


=U\  '  + 


5  —  y/E 
10 
Vs 


■A  tL{u* 


(3)  (2)  . 

U\  =U\  + 


A  tL{u^) 

A  tL{u'i>) 


5 

5  —  y/h 


10 


u?  =  (a)>n  +  b^AtL(un))  +  +  bY’AtLiu^))  +  (afjuf  +  b\z>AtL(uf) 

+  (a?lu(?)  +  b?)1AtL(u?)) 


,(0)„ 


do)  , 


(1) 


.(!)• 


l(2) 


.(2)> 
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a$un  +  b^AtL(un)J 
ah 2  4^  +  b^2AtL(u^ 

«S°3«n  +  &$A*Z(«n)) 

ai243^  +  b^3AtL(u^ 


a)>n  +  b^AtL(un)J 

ag-ug  +  bf^AtL(u^ 

a<2Au<2^  +  b2\AtL(u^ 
ag«n  +  b?lAtL(un)^ 
ag>«<3)  +  b^AtLiu? 
°235  m23)  +  b^lAtLiu:^ 


+  +  b^AtL(u^)j  +  +  b^2AtL(u[2)) 

))  +  {a2iu2)  +  b^AtL^)^ 

+  (agug  +  b^3AtL(u[L))^j  +  +  &gAiL(ug) 

))  +  +  tgAti^ug))  +  (agwg  +  b^3AtL(u(2 

+  (ai2Mi^  +  b{ilAtL(u[l)))  +  (agwg  +  b[2lAtL(u{2)) 


+  [a^u^  +  b^AtLiu^’) )  +  (agidr1  + 


,(2)„,(2)  .  j (2) 


ai,y  +  b{^AtL(un)J 

ag-ug  +  &g  AtZ^-ug 

«236  m23)  +  b^AtL^ 
ag  un  +  b^AtL(un)^ 

agwg  +  ^gAtl/^g 


(aggg  +  6^AtZ(ui1}) ) 
+  (ag^  +  b^lAtLiu^ 
+  (ag^g  +  b^AtLiu^ 
(ah6Mi1)  +  b^AtLiu^ 

+  (ag^g  +  &gAt4Mg 

+  (ggg  +  b^lAtL^ 


+  (aS42)  +  &b5A^(«iZ 


))  +  («225M22)  +  &225AtL(M22))) 


+  6gAhL(<0  J 

+  («gMg  +  b^2^AtL{%L2^ 


4242  +  b^AtL(u  g))  +  (agwg  +  ggA^g1 


+  («i26«i2)  +  b[2)iAtL(u(2))J 
))  +  («226m22)  +  b^AtL(ui2) 

))  +  («326M32)  +  &326AtL(M32) 

+  («129M12)  +  &SA^(MS2))) 

))  +  («229M22)  +  &229A^(M22) 

))  +  (429M32)  +  gWL(Mg 


+  («SM33)  +  &SAtL(M33))) 

with  the  coefficients  a  -  \,  and  gj,  given  by  (4.3),  (4.6),  (4.9),  (4.12),  (4.16),  (4.20), 
(4.34),  and 

4’i  =  0.2266,  42  =  0.2267,  42  =  0.2259,  42  =  0.1108, 

42  =  0.1602,  eng  =  0.2130,  42  =  0.1227,  42  =  0.0705, 

42  =  0.1409,  42  =  0.1290,  42  =  0.0911,  42  =  0.2797, 

(4.41)  42  =  0.6771,  42  =  0-1652,  42  =  0.0625,  /g1/  =  0.0923, 

=  0.2362,  /Zg  =  0.1739,  42  =  0.1772,  42  =  0.1370, 

42  =  0.2395,  42  =  0.0322,  /?g  =  0.0457,  /Zg  =  0.0724, 

$1  =  0.0567,  ctg  =  0.1469,  42  =  0-1162,  42  =  0-0609, 
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42 

=  0.0910, 

(3) 

«3,5 

=  0.2932, 

42 

=  0.0287, 

42  = 

0.0583, 

43 

=  0.0319, 

(3) 
a  4,1 

=  0.2821, 

42 

=  0.2265, 

42  = 

0.4054, 

42 

=  0.0745, 

d(3) 

Ml, 2 

=  0.1062, 

42 

=  0.0999, 

(3) 

74,1  = 

0.0146, 

(3) 

74,2  : 

=  0.0381, 

(3) 

74,3 

=  0.0203, 

0i  = 

0.8523, 

6-2  =  ] 

f  .0000, 

d3  =  0.8972,  d4 

=  1.0000 

The  CFL  coefficient  for  the  SSP  scheme  (4.40)-(4.41)  is  c  =  1.0319,  and  it  has  17 
evaluations  of  L  or  L.  It  is  therefore  slightly  more  efficient  than  the  scheme  (4.38)- 
(4.39).  We  could  also  reduce  the  number  of  L  or  L  evaluations  to  16,  however  the 
CFL  coefficient  reduces  to  c  =  0.6775,  which  is  not  impressive  at  all. 

Remark  4.2.  Our  analysis  is  based  on  the  choice  of  the  Chebyshev  Gauss-Lobatto 
nodes  as  the  subgrid  points  inside  the  interval  [tn,tn+1].  We  could  also  perform 
an  analysis  for  the  more  general  class  of  the  fourth  order  DC  scheme  in  which  the 
subgrid  points  are  placed  arbitrarily  subject  to  a  symmetry  constraint.  We  have 
performed  this  analysis  for  the  simple  case  of  Of,  =  0  for  all  k,  and  have  failed  to 
find  a  better  scheme  in  terms  of  the  SSP  property.  We  will  not  present  the  details 
here  to  save  space. 


5.  A  NUMERICAL  EXAMPLE 

In  this  section,  we  perform  a  numerical  study  to  assess  the  performance  of  the  DC 
time  discretizations,  coupled  with  the  fifth  order  weighted  essentially  non-oscillatory 
(WENO)  finite  difference  spatial  operator  with  a  Lax-Friedrichs  flux  splitting  [7],  to 
solve  the  following  Burgers  equation 

(5.1)  Mt+(y)  =0, 

with  the  initial  condition 

1  2 

(5.2)  u(x,  0)  =  -  +  -sin(7nr) 

o  o 

and  a  periodic  boundary  condition.  The  exact  solution  is  smooth  up  to  t  —  A,  then 
it  develops  a  moving  shock  which  interacts  with  the  rarefaction  waves.  We  use  the 
WENO  spatial  operator,  rather  than  the  TVD  spatial  operator,  since  the  former 
gives  better  accuracy  and  is  used  more  often  in  applications,  even  though  the  latter 
fits  better  the  theoretical  framework  of  this  paper,  being  rigorously  satisfying  the 
TVD  property  (1.5)  for  the  total  variation  semi-norm. 

In  Table  1  we  list  the  L 1  errors  and  the  numerical  orders  of  accuracy,  at  the  time 
t  =  0.2  when  the  solution  is  still  smooth.  We  compute  with  the  third  and  fourth 
order  SSP  DC  schemes  (3.22)-(3.23)  and  (4.38)-(4.39),  with  the  correct  incorporation 
of  the  operator  L,  and  with  the  original  third  and  fourth  order  DC  schemes  (3.1) 
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and  (4.1),  using  the  same  values  of  6k  but  without  using  the  operator  L.  For  this 
test  we  take  the  CFL  number  to  be  0.6,  that  is 

(5.3)  max  \u™\— —  =  0.6 

v  '  j  1  3 'Ax 

This  choice  is  based  on  the  heuristic  argument  that  the  spatial  WENO  operator  is  a 
high  order  generalization  of  the  second  order  generalized  MUSCL  scheme  [9],  which 
is  TVD  for  Erst  order  Euler  forward  time  discretization  under  the  CFL  condition 
rnaXj  =  0.5.  We  clearly  observe  in  Table  1  that  the  designed  order  of  accuracy 

is  achieved  or  exceeded.  The  other  SSP  schemes  in  sections  3  and  4  yield  similar 
errors.  We  do  not  present  their  results  in  order  to  save  space. 

Table  1.  L 1  errors  and  numerical  orders  of  accuracy.  Burgers  equa¬ 
tion  with  the  initial  condition  (5.2).  t  =  0.2 


Number  of 

DC3 

SSP  DC3 

DC4 

SSP  DC4 

cells 

L1  error 

order 

L1  error 

order 

L±  error 

order 

L1  error 

order 

20 

9.36E-4 

- 

1.27E-3 

- 

9.20E-4 

- 

1.35E-3 

- 

40 

4.78E-5 

4.29 

6.62E-5 

4.26 

4.27E-5 

4.43 

6.46E-5 

4.39 

80 

2.16E-6 

4.47 

2.82E-6 

4.55 

1.29E-6 

5.05 

2.11E-6 

4.94 

160 

1.81E-7 

3.58 

2.04E-7 

3.79 

5.38E-8 

4.58 

9.31E-8 

4.50 

320 

2.02E-8 

3.16 

2.07E-8 

3.30 

1.81E-9 

4.89 

3.21E-9 

4.86 

640 

2.48E-9 

3.03 

2.49E-9 

3.06 

4.40E-11 

5.36 

7.62E-11 

5.40 

When  t  =  0.6,  the  discontinuity  has  already  appeared.  We  plot,  in  Figure  1,  the 
solution  obtained  with  the  third  and  fourth  order  regular  DC  schemes  (3.1)  and 
(4.1),  and  SSP  DC  schemes  (3.22)-(3.23)  and  (4.38)-(4.39),  using  the  CFL  condition 
(5.3)  with  IV  =  40  equally  spaced  grid  points.  We  can  see  that  the  numerical 
solutions  are  indeed  non-oscillatory.  It  seems  that  for  this  test,  the  regular  DC 
schemes  without  using  the  operator  L  also  produce  non-oscillatory  results  for  the 
CFL  condition  (5.3). 

Finally,  we  would  like  to  numerically  assess  how  large  the  CFL  number  we  can 
take  and  still  maintain  stability.  We  compute  using  both  the  third  and  fourth  order 
SSP  DC  schemes  (3.22)-(3.23)  and  (4.38)-(4.39),  and  the  original  third  and  fourth 
order  DC  schemes  (3.1)  and  (4.1)  using  the  same  values  of  9k  but  without  using 
the  operator  L,  to  t  —  2,  with  N  =  160  equally  spaced  grid  points,  with  an  ever 
increasing  CFL  number.  In  Figure  2,  we  plot  the  L 1  errors  of  the  numerical  solution 
versus  the  CFL  number  for  the  third  order  (left)  and  fourth  order  (right)  schemes. 
We  observe  that  the  SSP  DC  schemes  are  indeed  stable  for  larger  CFL  numbers 
than  the  corresponding  regular  DC  schemes,  and  the  CFL  numbers  for  stability  are 
much  larger  than  the  theoretically  predicted  values  in  Theorems  3.1  and  4.1. 
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Figure  1.  Burgers  equation  with  the  initial  condition  (5.2).  t  =  0.6. 
N  =  40  equally  spaced  grid  points.  CFL  number  0.6.  Left:  third 
order  DC  schemes;  right:  fourth  order  DC  schemes.  Solid  line:  the 
exact  solution.  Circles:  SSP  DC  schemes.  Crosses:  regular  DC 
schemes. 


Figure  2.  Burgers  equation  with  the  initial  condition  (5.2).  t  =  2.0. 
N  =  160  equally  spaced  grid  points.  L 1  errors  (in  logarithmic  scale) 
versus  the  CFL  number.  Solid  lines:  regular  DC  schemes;  dashed 
lines:  SSP  DC  schemes.  Left:  third  order  schemes;  right:  fourth  order 
schemes. 


6.  Concluding  remarks 

We  have  studied  the  strong  stability  preserving  (SSP)  property  of  the  second,  third 
and  fourth  order  deferred  correction  (DC)  time  discretizations.  The  technique  of  the 
analysis  can  also  be  applied  in  principle  to  higher  order  DC  methods,  although  the 
algebra  becomes  more  complicated.  It  seems  that  the  DC  methods  do  not  have  as 
large  CFL  coefficients  as  the  Runge-Kutta  methods  for  the  SSP  property.  However, 
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since  the  DC  methods  can  be  easily  designed  for  arbitrary  high  order  accuracy, 
they  have  a  good  application  potential  and  the  analysis  for  their  SSP  property 
will  be  useful  for  their  application  to  solve  method  of  lines  schemes  for  hyperbolic 
conservation  laws. 
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